
*12345678901234567890123456789012345678901234567890123456789012345678901234567890
capture log close
clear all
set more off

*	************************************************************************
* 	File-Name: 	Replication_LossAnalysis.do
*	Log-file:	na
*	Date:  		01/09/2020
*	Author: 	Michaël Aklin, Johannes Urpelainen
*	Purpose:   	.do file to replicate:
*					Aklin, Michaël and Johannes Urpelainen
*					"Trials and Tribulations: Lost Energy Access Gains in Rural 
*					India" Energy for Sustainable Development.
*
*				Please see the paper for the sources of the data. All survey
*				data come from:
*					Aklin, Michaël, Chao-yo Cheng, Johannes Urpelainen, Karthik 
*					Ganesan, and Abhishek Jain. 2016. "Factors affecting 
*					household satisfaction with electricity supply in rural 
*					India" Nature Energy 1: 16170.
*	************************************************************************

*	************************************************************************
* 	A. Data
*	************************************************************************

* 	Modify the path to the data if needed!
use "DataReplication.dta.dta", clear

*	************************************************************************
* 	B. Analysis
*	************************************************************************

*	************************************************************************
* 	B.1 Summary stats
*	************************************************************************

*	
quiet tab m1_q23_edu, gen(educ_dummy) 
label var educ_dummy1 "No formal schooling"
label var educ_dummy2 "Up to 5th standard"
label var educ_dummy3 "Up to 10th standard"
label var educ_dummy4 "12th standard or diploma"
label var educ_dummy5 "Graduate and above"

quiet tab m1_q25_caste, gen(caste_dummy)
label var caste_dummy1 "Scheduled Caste"
label var caste_dummy2 "Scheduled Tribe"
label var caste_dummy3 "Other Backward Caste"
label var caste_dummy4 "General"
label var caste_dummy5 "Other caste"

quiet tab m1_q24_religion, gen(religion_dummy)
label var religion_dummy1 "Hindu"
label var religion_dummy2 "Muslim"
label var religion_dummy3 "Other religion"

*	Summary stats
preserve
keep if year == 2015
estpost sum felecloss logexpenses factor1 female_head educ_dummy* caste_dummy* religion_dummy* ///
	grid logpop c01_2001_dist_town if felecloss != ., d
esttab using "./Manuscript/Tables/summarystats.tex", ///
	cells("mean(label(Mean) fmt(2)) p50(label(Median) fmt(1)) sd(label(S.D.) fmt(1)) min(label(Min.) fmt(0)) max(label(Max) fmt(0)) count(label(Obs.) fmt(0))") ///
	nonumber label replace noobs
eststo clear
restore

*	Patterns
sum hhid if m2_q68_elec == 1 & l.m2_q68_elec == 1
sum hhid if m2_q68_elec == 1 & l.m2_q68_elec == 0
sum hhid if m2_q68_elec == 0 & l.m2_q68_elec == 1
sum hhid if m2_q68_elec == 0 & l.m2_q68_elec == 0

* Table 1
tab elec_cat // Columns 1-2
tab elecloss // Column 3

*	# of districts that had at least 1 hh lose electricity

preserve
collapse (max) elecloss, by(districts wave)
drop if wave==1
sum districts if elecloss==1 & wave==2
sum districts if elecloss==0 & wave==2
restore

tab elec_cat if m1_q8_state =="Bihar"
tab elec_cat if m1_q8_state =="Jharkhand"
tab elec_cat if m1_q8_state =="Madhya Pradesh"
tab elec_cat if m1_q8_state =="Odisha"
tab elec_cat if m1_q8_state =="Uttar Pradesh"
tab elec_cat if m1_q8_state =="West Bengal"

*	Change in satisfaction
tab dsatisfaction

*	Did cost explain loss?
replace m2_q781_elecunsatisfy_expensive = m2_q781_elecunsatisfy_expensive-1
replace m2_q812_ngrid_expensive_bill = m2_q812_ngrid_expensive_bill-1
replace m2_q811_ngrid_expensive_connect = m2_q811_ngrid_expensive_connect - 1
replace m2_q813_ngrid_unreliable = m2_q813_ngrid_unreliable-1
replace m2_q77_electrified_satisfaction = m2_q77_electrified_satisfaction-1
replace m2_q77_electrified_satisfaction = . if m2_q77_electrified_satisfaction<0

replace m2_q783_elecunsatisfy_quality = m2_q783_elecunsatisfy_quality-1
replace m2_q784_elecunsatisfy_repair = m2_q784_elecunsatisfy_repair-1
replace m2_q782_elecunsatisfy_unavail = m2_q782_elecunsatisfy_unavail-1

* Sources of insatisfaction
sum m2_q781_elecunsatisfy_expensive if year == 2015
sum m2_q781_elecunsatisfy_expensive if f.elecloss == 0 & elec2015==1
sum m2_q781_elecunsatisfy_expensive if f.elecloss == 1 & elec2015==1

sum m2_q782_elecunsatisfy_unavail if year == 2015
sum m2_q782_elecunsatisfy_unavail if f.elecloss == 0 & elec2015==1
sum m2_q782_elecunsatisfy_unavail if f.elecloss == 1 & elec2015==1

sum m2_q783_elecunsatisfy_quality if year == 2015
sum m2_q783_elecunsatisfy_quality if f.elecloss == 0 & elec2015==1
sum m2_q783_elecunsatisfy_quality if f.elecloss == 1 & elec2015==1

sum m2_q784_elecunsatisfy_repair if year == 2015
sum m2_q784_elecunsatisfy_repair if f.elecloss == 0 & elec2015==1
sum m2_q784_elecunsatisfy_repair if f.elecloss == 1 & elec2015==1

* Summary
sum m2_q68_elec if year == 2015, d
sum m2_q68_elec if year == 2018, d
sum m2_q69_elec_hrs if grid==1 & year == 2015, d
sum m2_q69_elec_hrs if grid==1 & year == 2018, d
sum m2_q69_elec_hrs if grid==0 & year == 2015, d
sum m2_q69_elec_hrs if grid==0 & year == 2018, d

*	************************************************************************
* 	C. Analysis
*	************************************************************************

*	************************************************************************
* 	C.1. Loss by type of supply
*	************************************************************************

tab m3_q83_light_main
sum f.elecloss if m3_q83_light_main == 1
sum f.elecloss if m3_q83_light_main == 2
sum f.elecloss if m3_q83_light_main == 3
sum f.elecloss if m3_q83_light_main == 4
sum f.elecloss if m3_q83_light_main == 5


*	************************************************************************
* 	C.2. Change in satisfaction
*	************************************************************************

eststo: reg dm3_q85_light_main_satisfy elecloss, cluster(m1_q11_village)
eststo: reg dm3_q85_light_main_satisfy elecloss i.districts, cluster(m1_q11_village)
eststo: reg dm3_q85_light_main_satisfy elecloss l.logexpenses l.factor1 l.debt l.female_head	///
	i.l.m1_q23_edu i.l.m1_q25_caste i.districts, cluster(m1_q11_village)
esttab using "./Manuscript/Tables/reg_sat.tex"	///
	,  booktabs label replace ///
 	nodepvars se(2) b(2) ///
	star(* 0.10 ** 0.05 *** 0.01) ///
	compress nogaps r2(2) sfmt(%9.2f %9.0f) ///
	nomtitles	///
	indicate("District FE = *districts*" "Expenses = *logexpenses*" "Wealth Index = *factor1*" "Debt = *debt*" "Female = *female_head*" "Education = *m1_q23_edu*" "Caste = *m1_q25_caste*", label($\checkmark$ ))	///
	noconstant nonotes  nobaselevels
eststo clear
	

*	************************************************************************
* 	C.2. Sources of loss
*	************************************************************************


eststo: reg f.elecloss logexpenses i.districts, cluster(m1_q11_village)
	estadd local range "$\{0;1\}$" , replace
	quiet sum felecloss if e(sample)==1
	estadd scalar mean = r(mean)
	estadd scalar sd = r(sd)
eststo: reg f.elecloss factor1 i.districts, cluster(m1_q11_village)
eststo: reg f.elecloss i.m1_q23_edu i.districts, cluster(m1_q11_village)
eststo: reg f.elecloss female_head i.districts, cluster(m1_q11_village)
eststo: reg f.elecloss ib4.m1_q25_caste i.districts, cluster(m1_q11_village)
eststo: reg f.elecloss i.m1_q24_religion i.districts, cluster(m1_q11_village)
eststo: reg f.elecloss grid i.districts, cluster(m1_q11_village)
eststo: reg f.elecloss logexpenses factor1 female_head	///
	i.m1_q23_edu ib4.m1_q25_caste i.m1_q24_religion grid logpop c01_2001_dist_town i.districts, cluster(m1_q11_village)
esttab using "./Manuscript/Tables/reg_loss.tex"	///
	,  ///
	booktabs label replace ///
 	nodepvars se(3) b(3) ///
	stats(N r2 mean sd, fmt(0 2 2 2 2) ///
	labels("N" "$ R^2$" "Mean of outcome" "Standard dev. of outcome")) ///
	star(* 0.10 ** 0.05 *** 0.01) ///
	compress nogaps r2(2) sfmt(%9.2f %9.0f) ///
	nomtitles	///
	indicate("District FE = *districts*", label($\checkmark$ ))	///
	noconstant nonotes  nobaselevels
eststo clear

eststo: logit f.elecloss logexpenses i.districts, cluster(m1_q11_village)
eststo: logit f.elecloss factor1 i.districts, cluster(m1_q11_village)
eststo: logit f.elecloss i.m1_q23_edu i.districts, cluster(m1_q11_village)
eststo: logit f.elecloss female_head i.districts, cluster(m1_q11_village)
eststo: logit f.elecloss ib4.m1_q25_caste i.districts, cluster(m1_q11_village)
eststo: logit f.elecloss i.m1_q24_religion i.districts, cluster(m1_q11_village)
eststo: logit f.elecloss grid i.districts, cluster(m1_q11_village)
eststo: logit f.elecloss logexpenses factor1 female_head	///
	i.m1_q23_edu ib4.m1_q25_caste i.m1_q24_religion grid logpop c01_2001_dist_town i.districts, cluster(m1_q11_village)
esttab using "./Manuscript/Tables/logit_loss.tex"	///
	,  booktabs label replace ///
 	nodepvars se(2) b(2) ///
	star(* 0.10 ** 0.05 *** 0.01) ///
	compress nogaps sfmt(%9.2f %9.0f) ///
	nomtitles	///
	indicate("District FE = *districts*", label($\checkmark$ ))	///
	noconstant nonotes nobaselevels eqlabels(none)
eststo clear


*	************************************************************************
* 	C.3. Loss of hours
*	************************************************************************

eststo: reg f.dhours logexpenses i.districts if elec2015==1, cluster(m1_q11_village)
	estadd local range "$[-24;24]$" , replace
	quiet sum fdhours if e(sample)==1
	estadd scalar mean = r(mean)
	estadd scalar sd = r(sd)
eststo: reg f.dhours factor1 i.districts if elec2015==1, cluster(m1_q11_village)
eststo: reg f.dhours i.m1_q23_edu i.districts if elec2015==1, cluster(m1_q11_village)
eststo: reg f.dhours female_head i.districts if elec2015==1, cluster(m1_q11_village)
eststo: reg f.dhours ib4.m1_q25_caste i.districts if elec2015==1, cluster(m1_q11_village)
eststo: reg f.dhours i.m1_q24_religion i.districts if elec2015==1, cluster(m1_q11_village)
eststo: reg f.dhours grid i.districts if elec2015==1, cluster(m1_q11_village)
eststo: reg f.dhours logexpenses factor1 female_head	///
	i.m1_q23_edu ib4.m1_q25_caste i.m1_q24_religion grid logpop c01_2001_dist_town i.districts if elec2015==1, cluster(m1_q11_village)
esttab using "./Manuscript/Tables/reg_hoursloss.tex"	///
	,  	///
	booktabs label replace ///
 	nodepvars se(2) b(2) ///
	stats(N r2 mean sd, fmt(0 2 2 2 2) ///
	labels("N" "$ R^2$" "Mean of outcome" "Standard dev. of outcome")) ///
	star(* 0.10 ** 0.05 *** 0.01) ///
	compress nogaps r2(2) sfmt(%9.2f %9.0f) ///
	nomtitles	///
	indicate("District FE = *districts*", label($\checkmark$ ))	///
	noconstant nonotes  nobaselevels
eststo clear


eststo: reg f.hoursloss_dummy logexpenses i.districts if elec2015==1, cluster(m1_q11_village)
eststo: reg f.hoursloss_dummy factor1 i.districts if elec2015==1, cluster(m1_q11_village)
eststo: reg f.hoursloss_dummy i.m1_q23_edu i.districts if elec2015==1, cluster(m1_q11_village)
eststo: reg f.hoursloss_dummy female_head i.districts if elec2015==1, cluster(m1_q11_village)
eststo: reg f.hoursloss_dummy ib4.m1_q25_caste i.districts if elec2015==1, cluster(m1_q11_village)
eststo: reg f.hoursloss_dummy i.m1_q24_religion i.districts if elec2015==1, cluster(m1_q11_village)
eststo: reg f.hoursloss_dummy grid i.districts if elec2015==1, cluster(m1_q11_village)
eststo: reg f.hoursloss_dummy logexpenses factor1 female_head	///
	i.m1_q23_edu ib4.m1_q25_caste i.m1_q24_religion grid logpop c01_2001_dist_town i.districts if elec2015==1, cluster(m1_q11_village)
esttab using "./Manuscript/Tables/reg_hourslossdummy.tex"	///
	,  booktabs label replace ///
 	nodepvars se(2) b(2) ///
	star(* 0.10 ** 0.05 *** 0.01) ///
	compress nogaps r2(2) sfmt(%9.2f %9.0f) ///
	nomtitles	///
	indicate("District FE = *districts*", label($\checkmark$ ))	///
	noconstant nonotes  nobaselevels
eststo clear

eststo: reg f.hoursloss_dummy_v2 logexpenses i.districts if elec2015==1, cluster(m1_q11_village)
eststo: reg f.hoursloss_dummy_v2 factor1 i.districts if elec2015==1, cluster(m1_q11_village)
eststo: reg f.hoursloss_dummy_v2 i.m1_q23_edu i.districts if elec2015==1, cluster(m1_q11_village)
eststo: reg f.hoursloss_dummy_v2 female_head i.districts if elec2015==1, cluster(m1_q11_village)
eststo: reg f.hoursloss_dummy_v2 ib4.m1_q25_caste i.districts if elec2015==1, cluster(m1_q11_village)
eststo: reg f.hoursloss_dummy_v2 i.m1_q24_religion i.districts if elec2015==1, cluster(m1_q11_village)
eststo: reg f.hoursloss_dummy_v2 grid i.districts if elec2015==1, cluster(m1_q11_village)
eststo: reg f.hoursloss_dummy_v2 logexpenses factor1 female_head	///
	i.m1_q23_edu ib4.m1_q25_caste i.m1_q24_religion grid logpop c01_2001_dist_town i.districts if elec2015==1, cluster(m1_q11_village)
esttab using "./Manuscript/Tables/reg_hourslossdummy_v2.tex"	///
	,  booktabs label replace ///
 	nodepvars se(2) b(2) ///
	star(* 0.10 ** 0.05 *** 0.01) ///
	compress nogaps r2(2) sfmt(%9.2f %9.0f) ///
	nomtitles	///
	indicate("District FE = *districts*", label($\checkmark$ ))	///
	noconstant nonotes  nobaselevels
eststo clear

*	************************************************************************
* 	C.4. Blackouts
*	************************************************************************


eststo: reg f.dblackouts logexpenses i.districts if elec2015==1, cluster(m1_q11_village)
eststo: reg f.dblackouts factor1 i.districts if elec2015==1, cluster(m1_q11_village)
eststo: reg f.dblackouts i.m1_q23_edu i.districts if elec2015==1, cluster(m1_q11_village)
eststo: reg f.dblackouts female_head i.districts if elec2015==1, cluster(m1_q11_village)
eststo: reg f.dblackouts ib4.m1_q25_caste i.districts if elec2015==1, cluster(m1_q11_village)
eststo: reg f.dblackouts i.m1_q24_religion i.districts if elec2015==1, cluster(m1_q11_village)
eststo: reg f.dblackouts grid i.districts if elec2015==1, cluster(m1_q11_village)
eststo: reg f.dblackouts logexpenses factor1 female_head	///
	i.m1_q23_edu ib4.m1_q25_caste i.m1_q24_religion grid logpop c01_2001_dist_town i.districts if elec2015==1, cluster(m1_q11_village)
esttab using "./Manuscript/Tables/reg_blackouts.tex"	///
	,  booktabs label replace ///
 	nodepvars se(2) b(2) ///
	star(* 0.10 ** 0.05 *** 0.01) ///
	compress nogaps r2(2) sfmt(%9.2f %9.0f) ///
	nomtitles	///
	indicate("District FE = *districts*", label($\checkmark$ ))	///
	noconstant nonotes  nobaselevels
eststo clear

eststo: reg f.dblackouts_dummy logexpenses i.districts if elec2015==1, cluster(m1_q11_village)
eststo: reg f.dblackouts_dummy factor1 i.districts if elec2015==1, cluster(m1_q11_village)
eststo: reg f.dblackouts_dummy i.m1_q23_edu i.districts if elec2015==1, cluster(m1_q11_village)
eststo: reg f.dblackouts_dummy female_head i.districts if elec2015==1, cluster(m1_q11_village)
eststo: reg f.dblackouts_dummy ib4.m1_q25_caste i.districts if elec2015==1, cluster(m1_q11_village)
eststo: reg f.dblackouts_dummy i.m1_q24_religion i.districts if elec2015==1, cluster(m1_q11_village)
eststo: reg f.dblackouts_dummy grid i.districts if elec2015==1, cluster(m1_q11_village)
eststo: reg f.dblackouts_dummy logexpenses factor1 female_head	///
	i.m1_q23_edu ib4.m1_q25_caste i.m1_q24_religion grid logpop c01_2001_dist_town i.districts if elec2015==1, cluster(m1_q11_village)
esttab using "./Manuscript/Tables/reg_blackoutsdummy.tex"	///
	,  booktabs label replace ///
 	nodepvars se(2) b(2) ///
	star(* 0.10 ** 0.05 *** 0.01) ///
	compress nogaps r2(2) sfmt(%9.2f %9.0f) ///
	nomtitles	///
	indicate("District FE = *districts*", label($\checkmark$ ))	///
	noconstant nonotes  nobaselevels
eststo clear


*	************************************************************************
* 	C.5. Type of off-grid
*	************************************************************************

eststo: reg f.elecloss ibn.light_num i.districts, cluster(m1_q11_village) noconstant
eststo: reg f.elecloss ib1.light_num i.districts, cluster(m1_q11_village)
	test i3.light_num +  i5.light_num == 0
	estadd scalar f = r(F)
	estadd scalar pf = r(p)
eststo: reg f.elecloss ib1.light_num logexpenses factor1 female_head	///
	i.m1_q23_edu ib4.m1_q25_caste i.m1_q24_religion logpop c01_2001_dist_town i.districts, cluster(m1_q11_village)
	test i3.light_num +  i5.light_num == 0
	estadd scalar f = r(F)
	estadd scalar pf = r(p)
esttab using "./Manuscript/Tables/reg_loss_offgrid.tex"	///
	,  ///
	booktabs label replace ///
 	nodepvars se(2) b(2) ///
	stats(N r2 f pf, fmt(0 2 2 2 4) ///
	labels("N" "$ R^2$" "F stat (Micro-grid $+$ SHS)" "p value for F stat")) ///
	star(* 0.10 ** 0.05 *** 0.01) ///
	compress nogaps r2(2) sfmt(%9.2f %9.0f) ///
	nomtitles	///
	indicate("Controls = logexpenses factor1 female_head *m1_q23_edu* *m1_q25_caste* *m1_q24_religion* logpop c01_2001_dist_town" ///
		"District FE = *districts*", label($\checkmark$ ))	///
	noconstant nonotes  nobaselevels
eststo clear


*	************************************************************************
* 	C.6. Quality of infrastructure
*	************************************************************************

label var m2_q781_elecunsatisfy_expensive "Electricity is too expensive"
label var m2_q782_elecunsatisfy_unavail "Electricity is not available when needed"
label var m2_q783_elecunsatisfy_quality "Electricity is of poor quality"
label var m2_q784_elecunsatisfy_repair "Maintenance and repair is poor"
label var m2_q77_electrified_satisfaction "Satisfaction with electricity"

eststo: reg f.elecloss m2_q77_electrified_satisfaction i.districts if grid==1, cluster(m1_q11_village)
	quiet sum m2_q77_electrified_satisfaction if e(sample)==1
	estadd scalar mean = r(mean)
eststo: reg f.elecloss m2_q781_elecunsatisfy_expensive i.districts if grid==1, cluster(m1_q11_village)
	quiet sum m2_q781_elecunsatisfy_expensive if e(sample)==1
	estadd scalar mean = r(mean)
eststo: reg f.elecloss m2_q782_elecunsatisfy_unavail i.districts if grid==1, cluster(m1_q11_village)
	quiet sum m2_q782_elecunsatisfy_unavail if e(sample)==1
	estadd scalar mean = r(mean)
eststo: reg f.elecloss m2_q783_elecunsatisfy_quality i.districts if grid==1, cluster(m1_q11_village)
	quiet sum m2_q783_elecunsatisfy_quality if e(sample)==1
	estadd scalar mean = r(mean)
eststo: reg f.elecloss m2_q784_elecunsatisfy_repair i.districts if grid==1, cluster(m1_q11_village)
	quiet sum m2_q784_elecunsatisfy_repair if e(sample)==1
	estadd scalar mean = r(mean)
eststo: reg f.elecloss m2_q781_elecunsatisfy_expensive m2_q782_elecunsatisfy_unavail m2_q783_elecunsatisfy_quality m2_q784_elecunsatisfy_repair i.districts if grid==1, cluster(m1_q11_village)
eststo: reg f.elecloss m2_q781_elecunsatisfy_expensive m2_q782_elecunsatisfy_unavail m2_q783_elecunsatisfy_quality m2_q784_elecunsatisfy_repair logexpenses factor1 female_head	///
	i.m1_q23_edu ib4.m1_q25_caste i.m1_q24_religion logpop c01_2001_dist_town i.districts if grid==1, cluster(m1_q11_village)
esttab using "./Manuscript/Tables/reg_loss_infra.tex"	///
	,  ///
	booktabs label replace ///
 	nodepvars se(2) b(2) ///
	stats(N r2 mean, fmt(0 2 2 2) ///
	labels("N" "$ R^2$" "Mean Indep. Var.")) ///
	star(* 0.10 ** 0.05 *** 0.01) ///
	compress nogaps r2(2) sfmt(%9.2f %9.0f) ///
	nomtitles	///
	mgroups("Grid Users" "Unsatisfied Grid Users", pattern(1 1 0 0 0 0 0 ) prefix(\multicolumn{@span}{c}{) suffix(}) span erepeat(\cmidrule(lr){@span}))		///
	indicate("Controls = logexpenses factor1 female_head *m1_q23_edu* *m1_q25_caste* *m1_q24_religion* logpop c01_2001_dist_town" ///
		"District FE = *districts*", label($\checkmark$ ))	///
	noconstant nonotes  nobaselevels
eststo clear


*	************************************************************************
* 	D. Graphs
*	************************************************************************


spagplot m2_q68_elec year	///
	, 	///
	id(districts)	///
 	xtitle("Year") ///
	note("")	///
	xlabel(2015 "2015" 2018 "2018")
graph export "./Manuscript/Figures/spagh_ElecAccess_district.pdf", replace as(pdf)


spagplot m3_q85_light_main_satisfy year	///
	, 	///
	id(districts)	///
 	xtitle("Year") ///
	note("") ///
	xlabel(2015 "2015" 2018 "2018")
graph export "./Manuscript/Figures/spagh_LightSatisfaction_district.pdf", replace as(pdf)


spagplot m2_q69_elec_hrs year	///
	, 	///
	id(districts)	///
 	xtitle("Year") ///
	note("") ///
	xlabel(2015 "2015" 2018 "2018")
graph export "./Manuscript/Figures/spagh_ElecHours_district.pdf", replace as(pdf)







